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ABSTRACT 

The treatment of dust attenuation is crucial in order to compare the predictions of 
galaxy formation models with multiwavelength observations. Most past studies have 
either used simple analytic prescriptions or else full radiative transfer (RT) calcu- 
lations. Here, we couple star formation histories and morphologies predicted by the 
semi-analytic galaxy formation model MORGANA with RT calculations from the spec- 
trophotometric and dust code GRASIL to create a library of galaxy SEDs from the 
UV/optical through the far Infrared, and compare the predictions of the RT calcula- 
tions with analytic prescriptions. We consider a low and high redshift sample, as well 
as an additional library constructed with empirical, non-cosmological star formation 
histories and simple (pure bulge or disc) morphologies. Based on these libraries, we 
derive fitting formulae for the effective dust optical depth as a function of galaxy phys- 
ical properties such as metallicity, gas mass, and radius. We show that such fitting 
formulae can predict the F-band optical depth with a scatter smaller than 0.4 dex 
for both the low and high redshift samples, but that there is a large galaxy-to-galaxy 
scatter in the shapes of attenuation curves, probably due to geometrical variations, 
which our simple recipe does not capture well. However, our new recipe provides a bet- 
ter approximation to the GRASIL results at optical wavelength than standard analytic 
prescriptions from the literature, particularly at high redshift. 
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1 INTRODUCTION 

Infrared (IR) observations of galaxies have demonstrated 
conclusively that dust is a fundamental and ubiquitous com- 
ponent of the interstellar medium (ISM) . Dust grains modify 
the chemical and physical conditions of the ISM by locking 
up a large fraction of heavy elements ejected by stars, by 
shielding molecular clouds (MCs) from dissociating radia- 
tion, thereby allowing them to cool and condense to the den- 
sities necessary to form stars, and by favoring the formation 
of H2 molecules themselves on their surfaces (see e.g. the re- 
views by Dorschner & Henning 1995; Draine 2003). There- 
fore dust grains are a fundamental ingredient in the very 
process of star formation. Also, dust plays a major role in 
shaping the spectral energy distribution (SED) of galaxies: 
dust grains absorb and efficiently scatter short wavelength 
(A 1/J.m) radiation. The extinction efficiency drops steeply 

for longer wavelength and dust can not survive at temper- 
atures £ 1000K: the absorbed energy is thermally emitted 



in the IR. Therefore dust modifies the intrinsic (pure stellar 
and/or non-thermal) SED of galaxies. This effect is partic- 
ularly strong in star-forming regions and starburst galaxies, 
where the youngest stars are deeply embedded within dense 
and optically thick dusty cocoons. 

The expected link between the intensity of the star 
formation activity and the amount of dust reprocessing in 
galaxies has been clearly observed both in the local and in 
the high redshift universe, revealing that much, and in many 
cases most, of the star formation activity is obscured in the 
Ultra-violet (UV) and optical, and can only be detected in 
the IR. For example, estimates based on observations with 
the IRAS satellite suggested that globally ~ 30% of the bolo- 
metric luminosity of nearby galaxies, mostly normal spirals, 
is reprocessed by dust in the IR (Soifer & Neugebauer 1991; 
Popescu & Tuffs 2002) . But IRAS also revealed the existence 
of a population of heavily obscured luminous and ultra- 
luminous IR galaxies (with Lir ~ 10 11 — 10 12 Lq, and with 
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L/ij £ 10 12 Lq, respectively, e.g. Sanders & Mirabel 1996), 
and provided a first hint of the strong evolution of these IR 
bright galaxies. This has been confirmed and quantified to 
high redshifts with surveys with ISO (z ~ 0.5—1, Elbaz et al. 
1999, 2002; Gruppioni et al. 2002; Dole et al. 2001), SCUBA 
(z ~ 2; Smail et al. 1997, 2002; Hughes et al. 1998; Chap- 
man et al. 2005) and Spitzer (z > 2; Le Floch et al. 2005; 

Babbedge et al. 2006). The optical to sub-millimetric cosmic 
background radiation provides a measure of the integrated 
star formation activity taking place at high-z. The energy 
density of the far-infrared (FIR) background measured by 
COBE (Puget et al. 1996; Hauser et al. 1998; Hauser & 
Dwek 2001) is comparable to that measured in the optical 
and near-IR (NIR) . All these observations indicate that the 
amount of energy emitted by dust over the history of the 
Universe is at least comparable to the energy emitted by 
stars that is able to escape galaxies and reach us, and that 
star formation activity and the consequent dust reprocessing 
is much more important at high redshift than locally. 

Therefore the presence and properties of dust in the 
ISM of galaxies, and its effect on the SEDs must be con- 
sidered with care in order to interpret observations to in- 
fer basic quantities such as star formation rates (SFR) and 
masses, as well as to compute reliable predictions (SEDs, 
luminosity functions, galaxy counts etc) from galaxy for- 
mation models. The various spectrophotometric codes (pe- 
GASE, Fioc & Rocca Volmerange, 1997; grasil, Silva et al., 
1998; STARBURST99, Leitherer et al., 1999; Bruzual & Char- 
lot, 2003) show reasonable agreement in predictions of the 
intrinsic, pure stellar SED of a galaxy, given its star forma- 
tion rate and metal enrichment history. However, the SED 
emerging from galaxies is the result of a complex interaction 
between the intrinsic properties of dust and the relative ge- 
ometry of the heating sources and the dust. The chemical 
composition, size distribution, and shape of the dust grains 
determine the absorption and scattering efficiencies of the 
dust mixture as a function of wavelength. These proper- 
ties are somewhat constrained only for our galaxy and the 
Magellanic Clouds through measurements of the extinction 
curve and the local diffuse dust emission. But it is also well 
known that these extinction curves have different shapes, a 
fact generally ascribed to different dust properties, and also 
within each galaxy the curves are spatially variable and de- 
pendent on the particular dusty environments sampled by 
the line of sight (e.g. Mathis et al. 1983; Rowan-Robinson 
1986; Cardelli et al. 1989; Mathis 1990; Fitzpatrick 1999, 
2007). A further issue is the degeneracy among different 
dust models that reproduce the average extinction curve and 
cirrus emission in the MW (e.g. Draine & Anderson 1985; 
Desert et al 1990; Dwek et al. 1997; Li & Draine 2001; Zubko 
et al. 2004). 

The relative geometry of stars and dust plays a major 
role in shaping the SED. This holds both for the attenuated 
starlight (e.g. Bruzual et al. 1988; Witt et al. 1992; Efs- 
tathiou & Rowan-Robinson 1995; Gordon et al. 1997; Fer- 
rara et al. 1999; Calzetti 2001) and for the dust emission 
spectrum. In fact, a fundamental property to be taken into 
account is that the relative geometry of stars and dust is age- 
and therefore wavelength-dependent, because the youngest 
stars, which dominate the UV luminosity, are also the most 
extinguished by the optically thick parent MCs (Silva et al. 



1998; Granato et al. 2000; Tuffs et al. 2004). In addition, 
detailed reproductions of the attenuation properties of disc 
galaxies seem to require an age-dependent extinction also 
in the diffuse medium for intermediate age with respect to 
older stars (Popescu et al. 2000; Panuzzo et al. 2007). All 
these factors then determine the temperature distribution 
of the dust grains, with each dusty environment and each 
dust grain having their own particular response to the ra- 
diation field, and therefore yielding the consequent shape of 
the emerging IR SED. 

In order to take full advantage of the wealth of multi- 
wavelength observations now available to us (UV/optical 
and IR), galaxy formation models must grapple not only 
with the modeling of stellar populations, but also with ab- 
sorption and re-emission by dust. This is equally true for 
both semi-analytic and numerical hydrodynamic simula- 
tions. In semi-analytic models (SAMs; e.g. White & Frenk 
1991), the evolution of the dark matter (DM) component is 
calculated directly using N-body methods or Monte Carlo 
techniques, while the evolution of the baryonic component 
is treated by simple recipes for the radiative cooling of gas, 
star formation, chemical evolution, feedback by supernovae 
and AGN, etc (see Baugh 2006 for a review). The SAM 
provides detailed information about the star formation and 
enrichment history of each galaxy, but typically only very 
limited information about the structural properties of the 
stars, gas, and dust (e.g., in general an effective radius can 
be computed for the disc and spheroid components) . 

Due to the complexities inherent to the treatment of 
the radiative effects of dust and the many unknowns con- 
nected to the dust properties, most SAMs have made use of 
simple empirical or phenomenological treatments. A widely 
adopted approach consists in computing the face-on dust 
optical depth at a reference wavelength (typically the V- 
band, Ty), and then computing the inclination dependence 
assuming that the stars and dust are uniformly mixed in a 
"slab" model, and that the wavelength dependence is given 
by a fixed "template" attenuation curve. Some modelers 
(e.g. Kauffmann et al. 1999; Somerville & Primack 1999; 
Nagashima et al. 2001; Mathis et al. 2002; De Lucia, Kauff- 
mann, & White 2004; Kang et al. 2005) compute ro,v using 
empirical relations between galaxy luminosity and dust op- 
tical depth (Wang & Heckman 1996). Others assume that 
t~o,v is proportional to the column density of dust in the 
disc, assuming that the dust mass is proportional to the 
metallicity (Guiderdoni & Rocca- Volmerange 1987; Lacey 
et al. 1993; Guiderdoni et al. 1998; Devriendt & Guiderdoni 
2000; Hatton et al. 2003; Blaizot et al. 2004; Kitzbichler 
& White 2007). De Lucia & Blaizot (2007) adopt the lat- 
ter prescription for obscuration due to the diffuse "cirrus" 
component, coupled with the approach suggested by Char- 
lot & Fall (2000) to account for the larger obscuration of 
the youngest stars by the molecular birth clouds. Cole et 
al. (2000) and Bell et al. (2003) instead couple their SAM 
with the Ferrara et al. (1999) library, which provides the 
net attenuation for smooth distributions of stars and dust 
as a function of wavelength and inclination angle based on 
Radiative Transfer (RT) calculations. 

In some works, the IR properties of galaxies in semi- 
analytic models are also computed. Using one of the above 
methods for computing the optical depth and attenuation 
curve, and hence the total amount of energy absorbed by 
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the dust, and assuming that all of this light is re-radiated 
in the IR, one can use observationally-calibrated templates 
describing the wavelength dependence of the dust emission 
(Guiderdoni et al. 1998; Devriendt, Guiderdoni, & Sadat 
1999; Devriendt & Guiderdoni 2000; Hatton et al. 2003; 
Blaizot et al. 2004), or modified Planck functions (e.g. Ka- 
viani et al. 2003) to compute IR luminosities. In a few cases, 
SAMs have instead been coupled with a self-consistent RT 
calculation for the full (UV to sub-mm) SED, e.g. Granato 
et al. (2000); Baugh et al. (2005), and Lacey et al. (2007, 
GALFORM); Fontanot et al. (2007, MORGANA); Granato et al. 
(2004) and Silva et al. (2005). 

Clearly there is a large range in both complexity and 
computation time represented in the treatments of dust in 
galaxy formation models in the literature. On the one hand, 
it is natural to be dubious that a simple empirical approach 
can adequately represent the intricacies of this very messy 
process. On the other hand, one of the advantages of SAMs 
is their flexibility and computational efficiency, and this ad- 
vantage is substantially undermined if one chooses to do full 
RT calculations for each galaxy. Moreover, given that SAMs 
are anyway based on a series of approximations, and that 
so many details of the relevant physics of galaxy formation 
(structural properties, star formation, chemical evolution, 
supernova and AGN feedback) are poorly understood and 
crudely modeled, feeding the results of SAMs into full RT 
may not be well motivated. 

In this paper, we make use of libraries of RT calculations 
based on SAM outputs as well as sets of non-cosmological, 
empirical star formation/enrichment histories to address the 
following questions: (1) How different are the results from 
the simple analytic dust models and the full RT calcula- 
tions? (2) How do the predictions from the RT calculations 
depend on the physical properties of the model galaxies? (3) 
Can we derive an improved analytic recipe, based on the RT 
calculations, that can be efficiently implemented in SAMs? 
The SAM that we use for this study is the morgana model 
(Monaco et al. 2007, Fontanot et al. 2006, 2007), and the 
RT libraries are constructed using the GRASIL code (Silva et 
al. 1998; Silva 1999). 

MORGANA produces good overall agreement with obser- 
vations of galaxy properties; any discrepancy is of little con- 
cern in this context, since in principle we could proceed by 
simply comparing SEDs resulting from different dust pre- 
scriptions applied to any possible realizations of the SAM. 
However, we wish to disentangle the dependence of the dust 
extinction predictions on galaxy physical properties from 
the correlations between these properties, and from distri- 
butions of these properties characteristic of a given cosmic 
epoch, as predicted by MORGANA or present in the real Uni- 
verse. Therefore we construct an "empirical" library, based 
on a classical "open box" chemical evolution code with sim- 
ple parameterized star formation histories, chosen to span 
the broadest range of possibilities. In addition, we consider 
both high-redshift and low-redshift catalogs extracted from 
the MORGANA outputs. In fact we may expect that locally 
calibrated simple dust prescriptions will give acceptable re- 
sults when coupled to a model characterized by a globally 
good agreement with local galaxy populations. Instead this 
may not be the case at high-z, where galaxies will typically 
have properties different from local ones. We focus on predic- 
tions of commonly used quantities in galaxy formation stud- 



ies, namely rest-frame magnitudes in the near-UV through 
Near-Infrared (SDSS ugri and 2MASS K-band). IR lumi- 
nosities from the mid-IR (3.5 /jm) to the sub-mm will be 
the focus of a companion paper (Fontanot et al., in prepa- 
ration) . 

The paper is organized as follows: in Section [2] we de- 
scribe the main features of the MORGANA and GRASIL models 
and the construction of the RT-based libraries that we will 
use for our comparison. In Section [3] we define some useful 
notation and summarize the analytic prescriptions for dust 
attenuation that we will evaluate. In Section [4] we derive 
fitting formulae between physical galaxy properties and the 
dust effective optical depth. In Section [SJ we present our 
analysis of the results of the RT calculations in terms of 
physical properties, and a comparison of the RT predictions 
at UV/optical/NIR wavelengths with several analytic pre- 
scriptions. We give our conclusions in Section[6] Throughout 
this paper we use magnitudes in the AB system (unless oth- 
erwise stated). 



2 GALAXY MODELS AND RT LIBRARIES 

The required inputs to the RT-based GRASIL model are the 
unextinguished stellar SED and the geometry and relative 
distribution of stars and dust. We assume that galaxies can 
be represented by a simple composite geometry of spheroid 
plus disc, and that both components are axisymmetric. The 
spheroid is represented by a King profile, and the disc by 
a radial and vertical exponential profile. The stellar SED is 
computed by convolving the distribution of stellar ages and 
metallicities arising from a given star formation and enrich- 
ment history with simple stellar population (SSP) models 
(Bressan, Granato, & Silva 1998; Bressan et al. 2002). The 
star formation histories, and hence stellar SEDs, are tracked 
separately for the bulge and disc components. We describe 
the details of how these quantities are obtained for each of 
our Libraries in the following sections. 

We make use of two different approaches to obtain en- 
sembles of these input parameters. The first approach is 
based on semi-analytic simulations with the MORGANA code 
(details given in H2.1l below) within the ACDM cosmological 
context. We extract catalogs at both low and high redshift 
from MORGANA. We refer to the libraries that are obtained 
from the MORGANA outputs as the Morgana Library (ML). 
In addition, we create a library in which the star forma- 
tion/enrichment histories are predicted by MORGANA as be- 
fore, but instead of having a composite disc+spheroid mor- 
phology, disc-dominated galaxies are represented as a pure 
disc and bulge-dominated galaxies as a pure bulge. In this 
way we can attempt to disentangle the effects of geometry 
and star formation history. We refer to this as the Control 
Library (CL). 

For the second approach, we use a classical chemical 
evolution code (che_evo), which is not embedded in a cos- 
mological context (details are given in i]2.2[) . The goal is 
to remove the "prior" that is contained in the cosmological 
models, in the form of the predicted correlations between 
and distributions of physical parameters. We refer to this as 
EL [Empirical Library). 

All libraries assume a Salpeter (1955) IMF with mass 
range from 0.1 to 100 Mq. Table [1] summarizes the differ- 
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Sample Name 




Star Formation History 


Geometry 


Redshift interval 


Morgana Library - Disc-dominated 


- Low-z 


MORGANA 


Disc+Bulgc 


0.0 < 


z 


< 0.2 


Morgana Library - Bulge-dominated 


- Low-z 


MORGANA 


Disc+Bulge 


0.0 < 


z 


< 0.2 


Morgana Library - Disc-dominated - 


High-z 


MORGANA 


Disc+Bulgc 


2.0 < 


z 


< 3.0 


Morgana Library - Bulge-dominated 


- High-z 


MORGANA 


Disc+Bulge 


2.0 < 


z 


< 3.0 


Control Library - Disc-dominated - 


Low-z 


MORGANA 


Pure Disc 


0.0 < 


z 


< 0.2 


Control Library - Bulge-dominated 


- Low-z 


MORGANA 


Pure Bulge 


0.0 < 


z 


< 0.2 


Control Library - Disc-dominated - 


High-z 


MORGANA 


Pure Disc 


2.0 < 


z 


< 3.0 


Control Library - Bulge-dominated - 


High-z 


MORGANA 


Pure Bulge 


2.0 < 


z 


< 3.0 


Empirical Library - Disc-dominated 


CHE_EVO 


Pure Disc 








Empirical Library - Bulge-dominated 


CHE_EVO 


Pure Bulge 









Table 1. Summary of RT libraries used in our comparison. 



ent libraries considered in this paper and described in this 
Section. 

2.1 Semi-analytic galaxy formation model: 

MORGANA 

MORGANA (MOdel for the Rise of GAlaxies aNd Active nu- 
clei) is a semi-analytic model for the formation and evolution 
of galaxies and AGNs (Monaco, Fontanot & Taffoni 2007; 
Fontanot et al. 2006, 2007). The main characteristics and in- 
gredients it implements are the following: (i) a sophisticated 
treatment of the mass and energy flows between galactic 
phases (cold and hot gas, stars) and components (bulge, disc 
and halo) coupled with the multi-phase treatment of the in- 
terstellar and intra-cluster media; (ii) an improved modeling 
of cooling and infall (Viola et al. 2007); (iii) a multi-phase de- 
scription of star formation and feedback processes (Monaco 
et al. 2004); (iv) a self-consistent description of AGN ac- 
tivity and feedback (Fontanot et al. 2006); (v) the building 
of the diffuse stellar component in galaxy cluster from the 
stars scattered in galaxy mergers (Monaco et al. 2006). 

The MORGANA realization from which we extract the li- 
brary of star-formation histories considered here, the ML, is 
the same presented in Fontanot et al. (2007), and we refer 
to the that paper for more details. It shows good agree- 
ment with the local stellar mass function, cosmic star for- 
mation history, the evolution of the stellar mass density, 
the slope and normalization of the Tully-Fisher relation for 
spiral discs, the redshift distributions and luminosity func- 
tion evolution of J^-band selected samples, and the 850/im 
number counts. Despite the fact that this model is able to 
reproduce the overall assembly of massive galaxies, it does 
not reproduce the "downsizing" trend of galaxies: in partic- 
ular it overpredicts the number of bright galaxies at z < 1 
(Monaco et al. 2006) and the number of faint galaxies at 
z 1 (Fontana et al. 2006). Moreover the model still under- 
predicts the counts of the brightest 850/im galaxies at high-z. 
In this model we neglect other ingredients, such as quasar- 
triggered galactic winds (see Monaco & Fontanot 2005). We 
checked that including these ingredients does not change the 
conclusions presented in this paper. 

MORGANA gives predictions for the star formation histo- 
ries of disc and bulge components separately, and for galaxy 
sizes. Disc sizes are computed using the Mo, Mao & White 
(1998) approach. The spin parameter of the DM halo is 



randomly chosen from its known distribution, and it is as- 
sumed that the angular momentum is conserved. The pres- 
ence of a central bulge is taken into account in computing 
the disc size. Bulge sizes are computed (as, e.g., in Cole et 
al. 2000) by assuming that kinetic energy is conserved when 
the merger of two collisionless discs takes place. We classify 
the objects in the ML into two categories according to the 
stellar bulqe-to-total mass ratio BIT = Mbu } 9 ^, . Fol- 

J ' M bidge+ A1 disc 

lowing Magdwick et al. (2003) we assume bulge-dominated 
galaxies to correspond to B/T > 0.6, and disc-dominated to 
B/T s: 0.6. 

We compare the physical properties of our model galax- 
ies with the SDSS sample, and in particular with the specific 
star formation r s f r = SFR/M+ vs. stellar mass M* plot pre- 
sented in Brinchmann et al. (2004) . They find that massive 
galaxies show low r s / r values, whereas less massive galaxies 
populate a horizontal strip in the diagram, corresponding to 
higher r s f r values. Moreover, Madgwick et al. (2003) find a 
relation between the morphological type and r s f r : objects 
with log(r rs /) < —11 (log(r rs /) > —11) tend to be bulge- 
(disc-) dominated galaxies. We show the corresponding plot 
for the objects in the ML in Fig. [1] we mark the position 
of bulge- (disc-) dominated galaxies as red crosses (cyan cir- 
cles). It is evident from this plot that the massive objects 
in MORGANA favor a bulge-dominated morphology, whereas 
the low-mass tail is populated by a large fraction of disc- 
dominated objects, in qualitative agreement with the SDSS 
observations. However, many massive galaxies show resid- 
ual star formation activity, which is not seen in the SDSS 
data: this is likely due to the incomplete quenching of cool- 
ing flows at late times in MORGANA (see Monaco et al. 2007 
for a complete discussion of this point). Moreover the star 
formation activity of many disc-dominated galaxies is too 
low with respect to observations, in part likely due to the 
strangulation processes acting on cluster satellite galaxies. 
However, these residual discrepancies between observations 
and model are not a major issue for this paper. 

Many relevant properties of galaxies, such as the star 
formation rates and sizes, evolve with redshift. We consider 
two samples drawn from MORGANA, a "low-z" sample at z < 
0.20 and a "high-z" sample at 2.0 < z < 3.0. In the following 
we will refer to the libraries based on these two samples 
as the "low-z" and "high-z" libraries. For both libraries we 
define disc- and bulge-dominated sub-samples. Starting from 
these samples we also define the CL accordingly. 
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Figure 1. log(M*/Af0) versus log(r s f r /yr 1 ) for the MORGANA library (ML). Red crosses and cyan circles refer respectively to bulge- 
dominated and disc-dominated objects. Left and right panel refer to the low-z and high-z sample respectively. 



Parameter 


lower value 


upper value 


^■wind 


1 Gyr 


15 Gyr 


h 


3 Gyr 


10 Gyr 


*2 


h 


*1 +7 Gyr 


Afinfall 


W 7 M Q 


1O 14 M 


T infall 


5 Gyr 


10 Gyr 


/* 





1 


e* 


Gyr 


7 Gyr 


e sch — en 





10 



Table 2. Parameters for the generation of random star formation 
histories with CHE_EVO. 



2.2 Simple star formation history model: CHE_EVO 

In order to explore the correlation between dust attenuation 
and the physical properties of model galaxies, we need to 
define a different sort of model library. The recipes that de- 
termine galaxy evolution in MORGANA (star formation, feed- 
back etc.) are responsible for built-in correlations among the 
physical quantities that enter in the determination of dust 
attenuation (metallicity, mass of gas and/or stars, star for- 
mation rate, radius of the galaxy, bolometric luminosity). 
These recipes (and their induced correlations) differ from 
SAM to SAM. In order to make sure our results are indepen- 
dent of the particular SAM we consider, we define another li- 
brary, EL, using the chemical evolution code CHE_evo (Silva 
1999). 

This is a "classical open box" chemical evolution code 
to compute the star formation and metallicity evolution 
of galaxies, with the following characteristics, (i) The gas 
from which stars form has an infall rate given by an ex- 
ponential function of time with an e- folding time Tinfaii. 



The normalization of the gas infall rate is fixed by set- 
ting the total mass Mjnfali at a given time tinfaii- (ii) The 
SFR can be expressed as a function of the mass of g £tS 5 ELS 
an arbitrary function of time, or as a combination of both: 
SFR(t) = z/ sch M gas (t) feach + /(£; A/ ga s), where the first term 
is a Schmidt-type law with efficiency f sc h and exponent k ac h- 
The second term represents a superimposed mode for which 
we set the initial and final time ti and t% and the fraction of 
gas converted into stars /*. We choose two possibilities for 
this mode: an exponential form depending only on time with 
e- folding time e*, and a "Schmidt enhanced" mode, with ef- 
ficiency e sc h-cn and exponent a sc h- cn - (iii) It is possible to 
set a time £ W md to represent the development of a galactic 
wind after which all the gas is removed, and the SFR and 
infall are stopped. All the gas present after this age comes 
from stellar evolution, (iv) The metallicity evolution Z(t) is 
computed by accounting for the lifetimes and yields of stars 
as a function of their mass and metallicity, including also 
type la supernovae. The evolution of all the quantities is 
then computed to a final age tfi n . 

In order to define a library with a uniform coverage 
of the possible evolutionary histories, we randomly extract 
the values of most of the Che_evo parameters within large 
intervals, as reported in Table [2] For each combination of 
those parameters we consider t nn = 5, 10, 13.6 Gyr. We fix 
^mfaii = 13 Gyr (which sets the normalization of the infall 
rate), fc scn = 1, a sc h-on = 1, since their variation does not 
add further information. In order to break possible remain- 
ing correlations, we then renormalize separately the star- 
formation history and the evolution of the cold gas mass to 
final values extracted randomly from uniform distributions. 
This choice forces the final metallicity evolution to be in- 
dependent from the star-formation history and the cold gas 
content. The final result of this procedure is a uniform dis- 
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tribution of the models in the (SFR, M gas , Z) space. We 
then randomly assign simple morphologies (i.e. pure disc or 
pure bulge) to the objects in the EL. We also recompute the 
EL using the Madgwick et al. (2003) criterion for assigning 
morphologies based on their specific star formation rates: 
we assume that objects with log(r s / r ) < — 11 correspond to 
pure spheroidal galaxies, and objects with log(r s / r ) ^ — 11 
are considered to be pure discs. We have checked that our 
conclusions do not change between the two libraries, so wc 
will present results only for the first realization. 



2.3 Dust and Radiative Transfer model: GRASIL 

For each object in the morgana and che_evo libraries, we 
compute the corresponding UV to radio SED using the spec- 
trophotometric code GRASIL (Silva et al. 1998). Subsequent 
updates and improvements are described in full detail in 
Silva (1999); Granato et al. (2000); Bressan, Silva & Granato 
(2002); Panuzzo et al. (2003); Vega et al. (2005). We refer 
to these papers for more details. 

GRASIL solves the equation of RT, taking into account 
a state-of-the-art treatment of dust effects, and includes the 
following main features: (i) stars and dust are distributed 
in a bulge (King profile) + disc (radial and vertical expo- 
nential profiles) axisymmetric geometry; (ii) the clumping of 
both (young) stars and dust through a two-phase interstel- 
lar medium with dense giant star-forming molecular clouds 
embedded in a diffuse ("cirrus") phase are considered; (iii) 
the stars are assumed to be born within the optically thick 
MCs and to gradually escape from them on a time-scale i csc , 
this gives rise to the age- (wavelength-) dependent extinc- 
tion with the youngest and most luminous stars suffering 
larger extinction than older ones; (iv) the dust composition 
consists of graphite and silicate grains with a distribution of 
grain sizes, and Polycyclic Aromatic Hydrocarbons (PAH) 
molecules; (v) at each point within the galaxy and for each 
grain type the appropriate temperature T is computed (ei- 
ther the equilibrium T for big grains or a probability distri- 
bution for small grains and PAHs) ; (vi) the radiative transfer 
of starlight through dust is computed along the required line 
of sight yielding the emerging SED; (vii) the simple stellar 
population (SSP) library (Bressan, Granato, & Silva 1998; 
Bressan et al. 2002) includes the effect of the dusty envelopes 
around AGB stars, and the radio emission from synchrotron 
radiation and from ionized gas in HII regions. 

For the MORGANA+GRASIL library (ML or CL), the scale 
radii for stars and dust in the disc (r* d , r d ) and in the bulge 
(rl, rjj) are provided by MORGANA. The disc scale heights 
for stars and dust (h* d and hf) are set to 0.1 times the cor- 
responding scale radii. 

For the CHE_evo+GRASIL library (EL), we assume that 
the scale radius depends on stellar mass as: 

r = krar {w^Wj kpc (1) 

where a r = 0.2 kpc for r£, a r = 3.0 kpc for r* d and k r is 
a random number between 0.5 and 2. We assume the same 
scale radii for the stellar and dust components and set the 
disc scale-heights to 0.1 times the corresponding scale radii. 

We fixed the other parameters needed by GRASIL and 
not provided by the star formation models, (i) The es- 



cape time-scale of young stars for the parent MCs is set 
to t osc = 10 7 yr, an intermediate value between that found 
by Silva et al. (1998) to describe well the SED of spirals 
(~ a few Myr) and starbursts (~ a few 10 Myr), and of 
order the estimated destruction time scale of MCs by mas- 
sive stars. In Fontanot et al. (2007) this value was found to 
produce good agreement with the K-band luminosity func- 
tions and the 850 /urn counts, (ii) The gas mass predicted 
by MORGANA or CHE_EVO is subdivided between the dense 
and diffuse phases, assuming the fraction of gas in the star 
forming molecular clouds /mc is 0.5. The results are not 
very sensitive to this choice, (iii) The mass of dust is ob- 
tained by the gas mass and the dust to gas mass ratio #dust 
which is set to evolve linearly with the metallicity given by 
the galaxy model, <5dust = 0.45 Z. (iv) The optical depth of 
MCs depends on the ratio tmc oc 8 duBt MmcAmc; we se t 
the mass and radius of MCs to typical values for the MW, 
Mmc = 10 6 Mq and tmc = 16 pc. (iv) The dust grain size 
distribution and composition is chosen to match the mean 
MW extinction curve. 

GRASIL computes the resulting SEDs along different 
lines of sight, but in the following we use only angle-averaged 
SEDs, unless otherwise stated. We convolve the SEDs with 
the GALEX FUV, SDSS u, g, r, i and 2-MASS K s filters. 

3 ANALYTIC AND EMPIRICAL 

PRESCRIPTIONS FOR ATTENUATION 

3.1 Definitions and notation 

For clarity, we first define some useful notation and the 
meaning of several quantities that we will use throughout 
the remainder of the paper. 

The amount and wavelength dependence of the dust at- 
tenuation on the intrinsic (pure stellar) SED of a system is 
commonly defined through the ratio between the emerging 
(after interaction with dust) and intrinsic luminosity, typi- 
cally expressed in magnitudes: 

A x = m x - ml = -2.5 log L x /L° x (2) 

It is important to clearly define and to explicitly dis- 
tinguish between the extinction and the attenuation curve 
(e.g. Granato et al. 2000; Calzetti 2001). The term extinc- 
tion is commonly used to describe the wavelength depen- 
dence of the optical properties (absorption plus scattering) 
of the dust mixture. Indeed, in the simple geometrical con- 
figuration of a slab of dust between the observer and a point 
source, as is the case when directly measuring the extinc- 
tion from observations of background stars, since the dust 
emission along the line of sight (both true dust emission and 
scattering) is negligible within a point source image, the so- 
lution of the equation of RT provides a direct link between 
the ratio of the observed to intrinsic energy, and the dust 
optical depth: L\/L x = exp(— t\). Therefore in this case 
A\ — 1.086 t\, where t\ = Tx, a bs + T\ tSca contains the infor- 
mation on the dust properties along the line of sight. 

This simple situation of course does not apply to studies 
of external galaxies. In this case we refer to A\ as the at- 
tenuation curve (other commonly used terms with the same 
meaning are obscuration, effective extinction or absorption) , 
since in this case the solution of the RT equation does not 
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provide a simple nor a priori predictable relation for L\ jL\ . 
For convenience, it is common practice to define an "effec- 
tive" optical depth such that L\jL\ — exp(—Tx ), with 
rtf^ summarizing the complex interaction between stars and 



dust in galaxies. This is further discussed in £15.11 

In the following, we refer to these quantities, the attenu- 
ation and the corresponding effective optical depth obtained 
with GRASIL as A^ s , t^ s , and the color excess as 



E GS (X 1 - Aa) 



(3) 



We will compare the values of these quantities obtained from 
GRASIL with those that we obtain when we apply various 
analytic dust prescriptions to the pure stellar SED. 



3.2 Analytic prescriptions for attenuation 

Several simple and computationally efficient treatments of 
dust effects in galactic SEDs have been proposed. They 
adopt different approaches to the problem, making a direct 
comparison difficult. However they also share some basic 
ingredients : (1) an attenuation law, either empirical (e.g. 
Calzetti et al. 1994; 2000) or obtained by coupling an ex- 
tinction law with an idealized geometrical configuration (2) 
a face-on effective optical depth at a reference wavelength 
(usually V-band), (3) a geometrical model for the depen- 
dence of the attenuation on viewing angle. Below we describe 
some of the most commonly adopted prescriptions for these 
ingredients, focusing on those that we will compare with in 
the next section. 



3.3 Extinction and Attenuation Laws 

For our own and a few nearby galaxies, the extinction law 
of the dust has been directly measured from observations of 
background stars. The differences found between the shapes 
of the (average) extinction curves of the Galaxy, the Large 
Magellanic Cloud and the Small Magellanic Cloud below 
~ 2600 A are often ascribed to the different metalicities in 
these systems. The Galactic extinction curve for the diffuse 
medium is commonly adopted (e.g. Mathis et al. 1983; Fitz- 
patrick 1989). 

Calzetti et al. (1994, 2000) have analyzed the dust at- 
tenuation for a sample of UV-bright starbursts. The derived 
attenuation curve is characterized by a shallower UV slope 
than that of the MW extinction, and by the absence of the 
2175 A feature. A complementary method of defining an 
empirical attenuation law has been presented by Chariot & 
Fall (2000, hereafter CF00). They consider an attenuation 
curve that is a power-law function of wavelength, with a nor- 
malization that depends on the age of the stellar population: 



CF 



r v (A/5500A)-°- 7 
AtT V (A/5500A)-°' 7 



if t ^ tMC 

if t > t M c 



(4) 



where tv indicates the total V-band attenuation experienced 
by young stars within the birth MCs (due to the birth clouds 
themselves and the diffuse ISM), tuc is the timescale of 
destruction of the MC or for stars to migrate out of the birth 
clouds and fi defines the fraction of the total effective optical 
depth contributed by the ambient ISM. CF00 used a sample 
of nearby UV-selected starburst galaxies to test their model 
and constrain its parameters, by considering the ratio of FIR 



to UV luminosities, the UV spectral slope, and Balmer line 
ratios (see also Kong et al. 2004). The best-fit values found 
by CF00 are t MC ^ W 7 yr, fj, ~ 0.3, with t v sm = fj,Tv ~ 0.5 
and t v ic = (l-^ry ~ 1.0. 

In the spirit of the CF00 model, one can of course adopt 
a composite attenuation curve, e.g. in which the Milky Way 
extinction curve is used to describe the diffuse "cirrus" dust 
and the CF00 power-law attenuation law is used to describe 
the situation experienced by young stars in the dense birth- 
clouds. This is the approach taken by De Lucia & Blaizot 
(2007, DLB07 hereafter). 

3.4 Effective Optical Depth 

Perhaps the most fundamental quantity characterizing the 
extinction is the face-on effective optical depth in a given 
reference band, usually the V band (tv). As we discussed in 
the introduction, many works have made use of the empiri- 
cal relation between the intrinsic luminosity in the B-band 
and the optical depth (Wang & Heckman 1996). Also very 
common are variants of the model originally proposed by 
Guiderdoni & Rocca-Volmerange (1987, hereafter GRV87), 
which links the optical depth of the dust at a given wave- 
length t\ to the column density of metals in the cold gas 
and to ry' 



2.1 x 10 21 cm" 2 



(•>) 



The exponent s of the metallicity is wavelength dependent: 
s = 1.35 for A < 2000A, and s = 1.6 for A > 2000A, based 
on interpolations of the extinction curves of the MW and 
the Magellanic Clouds. In GRV87 the column density Nh is 
assumed to be simply proportional to the cold gas fraction: 



{N H ) = 6.8 x 10 2 



M, 



A'L + M E 



(6) 



More recent models, in which an estimate of the radius 
of the disc has been made available, have used an updated 
form in which the column density is proportional to the cold 
gas mass divided by the gas radius squared (e.g. Guiderdoni 
et al. 1998; Devriendt & Guiderdoni 2000; Cole et al. 2000). 
DLB07 assume: 

AT M gas _ 2 

Nh = - , ;cm (7) 

1.4m p 7r r 2 - 

where r is the exponential scale radius of the disc, as pre- 
dicted by their semi-analytic model. 

As discussed above, the extinction suffered by young 
stars enshrouded in dense birth clouds may be greater than 
that experienced by more evolved stars embedded in more 
diffuse cirrus. We will adopt the age-dependent extinction 
prescription proposed by DLB07. They assume that stars 
younger than tuc (which they take to be 10 7 yr), have an 
optical depth r^ IC = r{> irrus (1/fx-l), with the value of fj, 
randomly extracted from a Gaussian distribution with (/x) = 
0.3 and o~„ = 0.2, truncated at 0.1 below and 1.0 above. 



3.5 Geometry 

The adopted extinction or attenuation curves are then gen- 
erally combined with simple geometrical configurations hav- 
ing analytical expressions, in order to describe the depen- 
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dence of the attenuation on the viewing angle i and/or in- 
clude scattering. The analytical solutions for the radiative 
transfer in simplified geometries are generally taken from 
Natta & Panagia (1984), Osterbrock (1989), Lucy et al. 
(1991), and Varosi & Dwek (1999). 

The simplest and most commonly used assumption is 
the slab configuration, where stars and dust are homoge- 
neously distributed in an infinite plane layer: 



Aff = -2.5 log 



1 - exp(-T A / cos(i)) 
t a /cos(i) 



(8) 



The t\ in the expression is sometimes "corrected" for 
scattering by multiplying it by VI — <*>\, where uj\ is the 
albedo (Guiderdoni & Rocca-Volmerange, 1987). Other ge- 
ometries, such as the oblate ellipsoid, have also been consid- 
ered (e.g. Devriendt et al. 1999). 



4 RELATING DUST EXTINCTION TO 
PHYSICAL PARAMETERS 

The GRV87 model for the optical depth, and its offshoots, 
provided a first basic attempt to relate the dust optical 
depth to physical properties of galaxies such as metallicity 
and gas column density. As we have discussed, empirically it 
is also known that dust optical depth correlates with proper- 
ties such as bolometric or B-band luminosity. Based on RT 
calculations with the ray-tracing code SUNRISE (Jonsson 
2006) applied within hydrodynamical simulations (Cox et 
al., 2004), Jonsson et al. (2006) found a tight relationship 
between effective optical depth and combinations of quan- 
tities such as metallicity, bolometric luminosity, star forma- 
tion rate, stellar mass or baryonic mass. Motivated by these 
results, we analyze our libraries to attempt to extract a gen- 
eralized fit for the effective optical depth as a power-law 
function of these variables, e.g. 



tv = /(iboi, M bar , SFR, Z, M gas , M* , r gal ) 



(9) 

As we discussed in ij2.2l correlations like the one found 
by Jonsson et al. (2006) may be in part the result of built- 
in correlations in the models, e.g., the Schmidt law used 
in the hydrodynamic simulations relates star formation rate 
to gas density. Therefore, we initially explore the depen- 
dence of dust obscuration on galaxy properties in our em- 
pirical CHE_evo+GRASIL library, EL, where we have deliber- 
ately broken all such built-in correlations. We then explore 
whether the same correlations hold in the ML. 

As we describe in !}2 in the EL, Z, M gas and SFR 
are uncorrelated by construction. The remaining 4 physical 
quantities (Lboi, M t , Mb ar and -R ga i) follow from the star 
formation history and our assumed correlation between M* 
and Rgai. We then perform a fitting procedure to try to 
recover the "true" effective optical depth from the GRASIL 
library, Ty S . We computed Ty S from the corresponding face- 
on GRASIL SEDs. We assume power-law dependences for the 
values of the 7 physical quantities and we define general 
relations involving independent physical quantities. We de- 
termine the best-fit parameters in each relation using a \ 2 
minimization technique and we repeat the fitting procedure 
considering the disc-dominated and the bulge-dominated 
subsamples separately. We start from a relation involving 
a single physical quantity, and then consider an additional 



degree of freedom if the resulting reduced \ 2 ls lowered by 
a factor ~ 2. Our results show that we usually need to in- 
volve at least three physical quantities in order to reproduce 



Ty well; and in particular the best fit result is obtained by 
considering the following relation: 





F 




T V 



AL 



M-. 



lkpc 



(10) 



In the following we refer to the ry values predicted by this 
fitting formula as Ty, and we indicate in parenthesis the 
library we use to calibrate the best-fit parameters. 

The optical depth Ty S in the EL is very clearly related 
to the gas mass M gas , scale radius 7? ga i and Z. This is not 
surprising, given that the dust-to-gas ratio we have used in 
GRASIL is proportional to the metallicity (see W1.2>\ . How- 
ever the relation between Ty S and Ty(EL) shows a large 
spread for low values of Ty S ; we checked that this behav- 
ior is connected to low values of the gas surface density 
E gas = M gas /(47ri?g aI ). For this reason, when computing the 
best-fit parameters for the EL (see tab. [3} we restrict our 
analysis to objects with S gas > 1Mqpc~ 2 for the EL discs 
and E flas > IOOMqpc -2 for the EL bulges (the different val- 
ues are related to the different meaning of the scale radius in 
the two geometries). The best-fit relation between Ty S and 
Ty (EL) is shown in Fig. [21 lower panels. Yellow dots refer to 
the whole EL, whereas blue and red dots refer to the objects 
with a gas surface density higher than the threshold. 

We then consider the SEDs in the CL. Recall that the 
CL uses the star formation and enrichment histories pre- 
dicted by MORGANA, but with bulge-dominated systems rep- 
resented by a pure bulge and disc-dominated systems repre- 
sented by a pure disc. We see that Ty (EL) provides a good 
approximation for Ty S in this library, for the disc sample. 
On the other hand most of the objects with bulge geom- 
etry lie below our surface density threshold, and therefore 
Ty (EL) is not a good solution for this subsample: in partic- 
ular it provides a systematic underestimate of the Ty S . We 
then redetermine the best-fit parameters for eq. [10] (i.e. the 
relation involving M gas , R ga i and Z) on the combined (low-z 
+ high-z) CL and collect the results in tab. [3] In fig. [5] (left 
panel) we show the comparison between Ty S and Ty(CL): 
for the pure discs the fit is accurate both for the combined 
and for the two separate samples; for the pure bulges the 
quality of the fit is slightly different for the "low-z" and 
"high-z" samples. 

When we include the effect of complex geometry in the 
ML, both Ty (EL) and Ty (CL) do not correctly reproduce 
Ty S . The CL fit parameters still provide a good approxima- 
tion for the intrinsic optical depth for the bulge-dominated 
subsamples. However, the presence of a bulge perturbs the 
SED of the disc-dominated objects. The effect grows when 
massive bulges are considered, and it is larger if the bulge 
component is actively forming stars. These effects are more 
important for the high-z galaxies, where both conditions 
are frequent in MORGANA. Again we re-calibrate the best 
fit parameters for eq. [10] on the combined (low-z + high-z) 
ML, and we collect our results in tab. [3] We then compare 
Ty (ML) and Ty S in fig. [2] (right panel). The figure shows 
that we obtain a reasonable fit to the combined sample, 
with a scatter of 0.1 dex and 0.4 dex for the disc and bulge 
samples respectively. However several problems arise when 
we consider the separate samples. We first discuss the disc- 
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Figure 2. Comparison between the "true" optical depth from the RT libraries Ty S and the predictions of our best-fit formulae feci, I10t . 
The bottom row shows the results from the Empirical Library (EL). The upper three rows show the results from the combined (high-z 
plus low-z), and separate high-z and low-z samples (as labeled) for the Control Library (CL; left plot) and the Morgana Library (ML; 
right plot). Vertical columns within each plot show the bulge-dominated, disc-dominated, and combined (bulge and disc) samples, as 
labeled. 





Sample 


Log(r°) 


a 


P 


7 


Scatter (dex) 


t?(EL) 


Discs 


-4.84 


0.39 


0.50 


-0.99 


0.2 




Bulges 


-3.63 


0.43 


0.34 


-0.65 


0.2 


t?(CL) 


Discs 


-5.59 


0.27 


0.56 


-1.05 


0.1 




Bulges 


-6.47 


0.21 


0.65 


-1.03 


0.5 


r(,(ML) 


Disc-Dom. 


-5.50 


0.19 


0.55 


-0.98 


0.1 


Bulge-Dom. 


-5.77 


0.32 


0.58 


-0.79 


0.4 



Table 3. Best-fit parameters for eq llOl for the different subsamples. For the CL and ML libraries the fits are for the combined low-z and 
high-z samples. 



dominated subsamples, where we find an overestimate at 
the level of ~ 1 dex for the high-z sample. This is due to 
objects with a significant bulge component (in mass, star- 
formation activity or both). On the other hand, when ap- 
plied to the bulge-dominated subsamples, Ty(ML) provides 
a systematic overestimate of Ty s for the high-z sample and 
a systematic underestimate for the low-z sample, both at 
the level of 0.5 dex. 

The fitting formula we propose in eq. [10] is analogous 
to standard analytic prescriptions, for example the GRV87 
and DLB07 approaches. However, it is interesting that the 
best-fit values for the parameters 5, e and r\ differ from the 
equivalent values generally used in the literature. Not too 
surprisingly, the combination of best-fit parameters for M gas 
and _R ga i hints of a dependence of Ty S on the cold gas surface 
density (for all fits in tab. [3] 7 ~ 2/3). This is similar to 
eq. even if in the latter equation the dependence is linear 
with the cold gas surface density. The dependence of Ty S on 



Prescription Optical Depth 



TAU-GS 
TAU-FIT 
TAU-GRV 
TAU-DLB 



"true" optical depth from GRASIL (Ty S ) 
results of fit to GRASIL libraries (t^(ML)) 
Guiderdoni & Rocca-Volmerange (1987; GRV87) 
De Lucia k, Blaizot, (2007; DLB07) 



Table 4. Dust Prescriptions 



metallicity Z is completely different, and much weaker than 
in eq. [5] 
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Figure 3. Left panel: The V-band normalized attenuation curves obtained by MORGANA+GRASIL for the low-z and high-z samples (yellow 
shaded area and diagonal texture respectively), are shown compared with the MW extinction law (dotted red line, by Mathis et al. 1983), 
the Calzetti et al. (dashed blue line), the CFOO (dot-dashed green line), and the slab + with the MW extinction (continuous red line) 
attenuation curves. Right panel: the same quantities for Prescription TAU-GS (see sect. 15. 2t . 



5 COMPARISONS AT OPTICAL 
WAVELENGTHS 

In this section, we compare the results of the GRASIL libraries 
with several approximations. In all of these approximations, 
we assume that the wavelength and geometry dependent be- 
havior of the dust extinction can be captured by using an 
expression that relates the effective optical depth tv to phys- 
ical properties of the galaxy, in combination with a "slab" 
model for the geometry and the DLB07 composite MW + 
CFOO attenuation curve (see following DLB07, we select 
/x randomly from a Gaussian distribution) . We assume a sin- 
gle dust slab to describe the attenuation of both disc and 
bulge. To assure proper comparison with the GRASIL SEDs, 
we also average slab model results over inclination. We will 
refer to the four prescriptions that we will test: 

• Prescription TAU-GS: we use the actual value of the 
effective optical depth Ty s from the RT library. 

• Prescription TAU-FIT: we use the fitting formula 
Ty (ML) derived in §H 

• Prescription TAU-GRV: we use the GRV87 formula 
(Eqn. [SJ for TV- Recall that in this approximation, the col- 
umn density of the gas is simply proportional to the gas 
fraction (no scaling with galaxy radius). 

• Prescription TAU-DLB: we use the DLB07 formula 
(Eqn. [7} for Ty. Here, the optical depth is proportional to 
the surface density of the gas (scales as l/r| al ). 

These prescriptions are summarized in Table [4] 

We have extensively tested many different possible 
combinations of ingredients, for example, substituting the 
Calzetti or CFOO attenuation law for the Galactic extinc- 
tion law, or substituting the oblate ellipsoid geometry for 



the slab model. We find a similar level of agreement with 
the Calzetti law, but worse agreement with the CFOO atten- 
uation law. The oblate ellipsoid geometry produces similar 
results to the slab model. 



5.1 Attenuation curves 

In Fig. [3] (left panel) we show the V-band normalized atten- 
uation curves Af s for the GRASIL SEDs in the ML. The yel- 
low shaded area is for the low-z sample, the diagonal texture 
for the high-z sample (l-cr confidence regions). It is worth 
noting that the galaxy-to-galaxy scatter shown in fig. [3] is 
an underestimate, as we do not consider the distribution 
of inclination angles (but only angle- averaged SEDs) and 
we keep fixed both t asc and the intrinsic properties of dust 
grains. For comparison we also display the average MW ex- 
tinction curve (by Mathis et al. 1983, dotted red line), and 
the attenuation curves of Calzetti et al. (2000) (dashed blue 
line), CFOO (dot-dashed green line), and the slab + MW ex- 
tinction + scattering (continuous red line). All curves have 
been normalized to Tv = 1 in order to highlight the different 
dependencies on the wavelength. 

The comparison with the Calzetti et al curve is par- 
ticularly interesting, since its shape has been interpreted in 
several ways. For instance, Gordon, Calzetti, & Witt (1997) 
ascribed the observed shape to the presence of dust lack- 
ing the 2175 A feature in the extinction curve, i.e. intrin- 
sic dust differences. This interpretation is probably a conse- 
quence of adopting a model that accounts for the clumping 
of dust but not of stars, and with a spatial distribution of 
stars independent of stellar age. Granato et al. (2000) repro- 
duced this featureless and shallow attenuation as a result of 
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the complex and wavelength-dependent geometry where the 
UV-emitting stars are heavily embedded inside MCs, while 
older stars, mainly emitting in the optical and NIR, suffer a 
smaller effect from the diffuse medium. Similar conclusions 
on the strong effect of complex geometries in the resulting 
attenuation have been found in the analysis by e.g. Panuzzo 
et al. (2003, 2007), Pierini et al. (2004), Tuffs et al. (2004), 
Inoue (2005), and Inoue et al. (2006) . While one can be con- 
fident that intrinsic dust differences are effectively present 
in different galaxies, as well as in different environments, it 
also appears that the complex geometry in galaxies can have 
a major and sometimes dominant role. 

Since, as described in 32.31 in the MORGANA+GRASIL 
models the dust properties (size distribution and composi- 
tion) are fixed so as to reproduce the average MW extinc- 
tion curve, any difference between that and the attenuation 
curves we obtain, must be ascribed to the particular ge- 
ometrical configuration of each model galaxy. As we said, 
the geometry has to be understood in a broad sense, i.e. 
not only the stars and dust distributions in a bulge and/or 
in a disc, but also the clumping of stars and dust and the 
age-, therefore wavelength-dependence of the amount of at- 
tenuation suffered by stars. We stress that the attenuation 
curve we obtain for each model would be exactly equal to 
the adopted extinction curve if we were to assume that all 
the dust is in a foreground screen in front of the stars and 
no scattered light reached the observer. Therefore the very 
large spread we get, as evident from the figure, has to be 
ascribed to the accounting of the star formation histories, 
the distribution of stars and dust with the appropriate radii 
and the connection between stellar populations of different 
ages with different dusty environments. 

We also show in fig. [3] (right panel) the correspond- 
ing mean attenuation as a function of wavelength for Pre- 
scription TAU-GS; i.e., using the composite attenuation law 
suggested by DLB07 based on the CF00 model. Encour- 
agingly, the behavior is qualitatively similar to that of the 
the intrinsic GRASIL mean attenuation curve. The major dis- 
crepancies arise at very short wavelengths (in both samples) 
where GRASIL predicts larger attenuations, and in the near- 
infrared region (in the high-z sample) , where GRASIL predicts 
a smaller scatter and on average larger attenuations. We test 
that the scatter in the attenuation laws at short wavelengths 
is reduced if we assume a fixed value for the fi parameter 
(i.e. fi — 0.3) in Prescription TAU-GS, instead of using a 
Gaussian distribution. 

5.2 Attenuation and color excess 

In this sub-section we show the results we obtain for atten- 
uation and reddening by coupling MORGANA either with the 
full GRASIL RT calculation, or by applying different prescrip- 
tions for attenuation to the pure stellar SED (summarized in 
Table f3| . We illustrate our comparison by showing attenua- 
tions and color excesses as a function of M* and specific star 
formation rate SFR/M*. We show the results for the FUV 
GALEX filter, the g and r SDSS filters and for the A"-band 
bandpass, but similar results holds if we consider the full set 
of SDSS bandpasses plus near-infrared photometry. We also 
consider bulge- and disc-dominated subsamples separately. 
We consider restframe magnitudes both for the low-z and 
high-z samples. This choice allow us to test the robustness 



of different prescription for dust attenuation, when changing 
the properties of the underlying galactic populations. 

In Fig. |4]we compare E GS (g - r), E GS (r - K), Apuv 
and A GS to the corresponding quantities obtained when r GS 
is combined with a composite attenuation law like DLB07 
(Prescription TAU-GS; shaded areas and diagonal texture 
represent the 1—a confidence regions). Clearly, this case pro- 
vides a reference point for whole analytic approach — in 
principle, none of our other prescriptions can possibly do 
better at reproducing the full GRASIL results than this one, 
where the effective optical depth is measured directly from 
the RT libraries. 

Indeed, we see that this model provides a satisfactory 
representation of the RT results at optical wavelengths. 
Mean attenuations are slightly underpredicted, but the 
trends as a function of M* and SSFR are correctly re- 
covered. However the agreement worsens at shorter wave- 
lengths, with a significant underprediction of the level of 
attenuation. This result is expected: at FUV wavelenghts 
the attenuation is extremely sensitive to the details of the 
dust distribution and composition, and it is very difficult to 
describe it with a simple analytic formulation. This result 
can also be understood in terms of the different distribu- 
tions of attenuation that we have already seen in Fig. [3] 
The same prescription is able to reproduce A G fj V and 
A GS for the high-z sample, but the agreement is poor for 
both color excesses, particularly in the disc-dominated sub- 
sample. Again, this reflects the disagreement in the attenu- 
ation curves. 

We compare the same quantities for our TAU-FIT pre- 
scription in fig. [5] The agreement is again satisfactory, and 
more importantly where there are discrepancies, they are 
similar to those seen in Prescription TAU-GS. This indi- 
cates that the main source of remaining error is not our 
fitting results for ry but rather the simplified approach we 
have used to estimate the attenuation curve. 

We show in fig. [S] the same comparison with Prescrip- 
tion TAU-GRV. In the low-z sample (left panel) , the ranges 
covered by the model and the prescription are in reasonable 
agreement, despite a slight systematic overestimate of atten- 
uation and color excess in both the disc- and in the bulge- 
dominated sample. The dependence of attenuation on M* is 
at least approximately reproduced. The dependence of at- 
tenuation and color excess on M* are only marginally repro- 
duced, with a larger discrepancy for the bulge-dominated ob- 
jects. Moreover, the distributions of attenuations and color 
excesses around the mean values show a smaller (larger) 
scatter with respect to GRASIL predictions for disc- (bulge- 
dominated objects. The effect of considering a simple ge- 
ometry is similar to the low-z case (higher attenuations and 
color excesses for the pure discs and vice- versa for the pure 
bulges). 

Finally, we show in ng.[7]the comparison with Prescrip- 
tion TAU-DLB. Following the original work, for each object 
we compute r^/i starting from the disc scale-radius (even 
if the galaxy is bulge-dominated). Prescription TAU-DLB 
is able to reproduce both GRASIL attenuations and color ex- 
cesses and their dependence on M* in the low-z sample, with 
the notable exception of the far-UV attenuations. Moreover, 
the scatter is always larger than the corresponding quantity 
in GRASIL. Quite interestingly a very similar behavior is ob- 
tained for the CL, with only small differences with respect 
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Figure 4. Comparison between the color excesses E GS (g — r) and E GS (r — K), and the attenuations A^^y and A^ s (shaded areas)as 
a function of stellar mass M t , (left panel) and of specific star formation rate (right panel), with the predictions of Prescription TAU-GS 
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labeled on the figure. 
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Figure 5. Same as in fig. [4] The vertical textured area shows the predictions of Prescription TAU-FIT. 



to the ML. We then conclude that Prescription TAU-DLB 
is less sensitive to the composite vs simple geometry than 
Prescription TAU-GRV. On the other hand, when we com- 
pare Prescription TAU-DLB to the high-z sample we obtain 
both a steeper dependence of attenuation from and a 
larger scatter with respect to GRASIL. 



Both Prescriptions TAU-GRV and TAU-DLB provide 
acceptable agreement with the GRASIL results for the low-z 
samples, but they do not do well at high redshift. The satis- 
factory agreement at low-z is not surprising, given the fact 
that the original prescriptions are calibrated to reproduce 
the properties of galaxies in the local Universe. The origin 
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Figure 6. Same as in fig. [4] The vertical textured area shows the predictions of Prescription TAU-GRV. 
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of the high-z discrepancy lies in the different distribution 
of cold gas fractions and surface densities in the two sam- 
ples (the dependence on Z and M gas being the same). The 
DLB07 prescription for ry is proportional to the gas surface 
density; we demonstrate that this assumption is broadly con- 
sistent with the RT results, however our analysis suggests a 
much weaker dependence of tv on surface density (to the 
square-root power). We therefore expect the differences to 



grow at high-z, where gas surface densities are considerably 
higher. 

From fig. [S] to [TT] we summarize the results for all four 
prescriptions in terms of the residuals of the difference be- 
tween the mean attenuation predicted by the prescription 
and the intrinsic attenuation in GRASIL. From these figures, 
we can clearly see that prescription TAU-GS always has 
the smallest residuals, as expected. Prescription TAU-FIT is 
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Figure 8. Residuals in the mean attenuations as a function of stellar mass for our four prescriptions with respect to GRASIL predictions. 
Dotted, dashed, dot-dashed, long-short dashed lines refer to prescription TAU-GS, TAU-FIT, TAU-GRV and TAU-DLB respectively. 



seen to be quite successful in that in most cases, the resid- 
uals are similar to those of prescription TAU-GS and the 
trends are similar. Prescriptions TAU-GRV and TAU-DLB 
show a similar level of agreement overall, but often show dif- 
ferent trends with galaxy properties such as stellar mass. All 
prescriptions seem to have particular difficulty reproducing 
the GRASIL results for high-z bulge-dominated galaxies. 

Finally, we assess the impact of the different dust pre- 
scriptions on the predictions of the statistical properties of 
galactic samples. In particular, in fig. \T%\ we focus on the 
resulting luminosity functions, which are commonly used 
as a baseline prediction for semi-analytic models. The thin 
solid line represents the unextinguished luminosity function 
at different wavelengths for the low-z and high-z samples, 
whereas the thick solid line shows the luminosity function 
of the samples as predicted by GRASIL. The other lines repre- 
sent Prescription TAU-GS, TAU-FIT, TAU-GRV and TAU- 
DLB (the dotted, dashed, dot-dashed and long-short dashed 



line respectively). This comparison is quite encouraging, in 
that it indicates that, while it is extremely important to 
correct for the effects of dust, particularly at high redshift 
and in shorter wavelength bands, even a simple analytic 
approach gives quite good results, when dealing with inte- 
grated quantities, such as the LFs, in which galaxy to galaxy 
differences are smoothed out, particularly in the optical-NIR 
region characterized by relatively regular and featureless 
SEDs. However, the spread in the predicted LF increases 
at UV wavelengths. As we have seen, the effect of dust on 
galaxy colors is not reproduced as well by the simple analytic 
prescriptions. 



5.3 Inclination dependence of attenuation 

The geometrical dependence of the predicted SEDs man- 
ifests itself also in the dependence of attenuation on the 
viewing angle. In many papers, statistical estimates of the 
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Figure 9. Residuals in the mean attenuations as a function of specific star formation rate for our four prescriptions with respect to 
GRASIL predictions. Lines as in fig. \E\ 



dust opacity in galaxy discs have been worked out, through 
the analysis of the inclination dependence of apparent mag- 
nitudes and colors of large samples of spiral galaxies (e.g. 
de Vaucouleurs et al. 1991; Giovanelli et al. 1995; Moriondo, 
Giovanelli & Haynes, 1998; Tully et al. 1998; Graham 2001; 
Masters, Giovanelli & Haynes, 2003). Extensive work on 
the angle dependence predicted by RT models (chiefly for 
disc galaxies) has also been presented in the literature (e.g. 
Kylafis & Bahcall 1987; Byun, Freeman & Kylafis 1994; 
Bianchi, Ferrara & Giovanardi 1996; Kuchinski et al. 1998; 
Tuffs et al. 2004; Pierini et al. 2004; Rocha et al. 2008). 
The interpretation of observations to infer opacities, scale- 
lengths, luminosities etc. are non-trivial, because of the ra- 
dial gradients in the dust distribution and metallicity thus 
in the dust opacity, the different scale-heights of stars and 
dust, the bulge to disc ratio, and the inclination effects (e.g 
Xilouris et al. 1999; Boissier et al. 2004; Popescu et al. 2005; 
Mollenhoff, Popescu & Tuffs 2006; Driver et al. 2007). 



The statistical data for the attenuation-inclination de- 
pendence may therefore provide a further test for SAM pre- 
dictions (see e.g. Granato et al. 2000), if the spectral prop- 
erties computed for the SAM are reliable, given the com- 
plexities described above. In the previous Sections we have 
considered the results we obtain for angle-averaged SEDs. In 
Fig. [T3]we show the dust net attenuation (difference edge-on 
to face-on) predicted using GRASIL and Prescription TAU- 
DLB, which we recall adopts an infinite plane parallel slab 
(Eq. [BJ . It is evident that the predicted net attenuations are 
very different in the two cases. We do not want to discuss 
here the agreement of the result with respect to data, we 
only wish to point out with one example the clearly dif- 
ferent behavior both in terms of trend and spread. Part of 
the difference is surely due to the fact that with a slab the 
net attenuation is independent from the bulge to disc ra- 
tio given by the SAM. Instead, as discussed by Tuffs et al. 
(2004), Pierini et al. (2004) and Driver et al. (2007), the 
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Figure 10. Residuals in the variance in the mean attenuations as a function of stellar mass for our four prescriptions with respect to 
GRASIL predictions. Lines as in fig. \E\ 



attenuation suffered by a dust-less bulge seen through the 
disc can be very pronounced as a function of the viewing 
angle, depending on the scale radii of stars and dust. In 
MORGANA most bulges host relevant star formation activity 
within dusty MCs, therefore in these cases the total effec- 
tive attenuation of the galaxy is dominated by (spherically 
symmetric) bulges rather than by discs, and this results in 
lower values and spread of the net attenuation as compared 
to a slab distribution. 



6 SUMMARY 

In this paper, we attempt to better understand and quan- 
tify the results of detailed radiative transfer and dust mod- 
els (grasil), with the goal of assessing and improving the 
treatment of dust in semi-analytic models of galaxy for- 
mation. We build several different libraries of star forma- 
tion/chemical enrichment histories, which we couple with 



the GRASIL RT code. One library is based on the MORGANA 
semi-analytic model, with galaxies extracted from two dif- 
ferent redshift intervals [z < 0.20 and 2 < z < 3). As a first 
choice, we keep the composite geometry (bulge+disc) pre- 
dicted by the SAM; we additionally define a Control Library 
by considering only the dominant galactic component (bulge 
or disc). We then build a third library using the chemical 
evolution code CHE_EVO and requiring that cold gas mass, 
metallicity and star formation rate are both independent 
and uniformly distributed in the (SFR, M gas , Z) space. For 
this Empirical library we consider a simpler geometry (pure 
bulge or pure disc). We then interface the star formation 
histories with the spectrophotometric code GRASIL to ob- 
tain the corresponding synthetic extinguished and unextin- 
guished SEDs. 

We compute the intrinsic optical depth in the l/-band 
Ty S for each object, and we combine this estimate with dif- 
ferent prescriptions for the attenuation law and for the rel- 
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Figure 11. Residuals in the variance in the mean attenuations as a function of specific star formation rate for our four prescriptions 
with respect to GRASIL predictions. Lines as in fig. [8] 



ative geometry of stars and dust. We find that the GRASIL 
results are best reproduced by a composite attenuation law 
as proposed by De Lucia & Blaizot (2007), which consists of 
a power-law in wavelength for the young stellar population, 
as suggested by Chariot & Fall (2000), and a Galactic ex- 
tinction curve plus a slab model for older stars. However, we 
notice that even a simpler model with a Galactic extinction 
law plus a slab model is already a good approximation for 
the GRASIL results, if the "true" value of r« s is known. 

We then try to estimate the dependence of Ty S on phys- 
ical properties of the model galaxies in the Empirical Li- 
brary. We find that Ty S depends mostly on the gas metal- 
licity, on the cold gas mass and on the radius of the sys- 
tem. We provide fitting formulae that are able to predict 
Ty S over a large dynamical range and with reasonable scat- 
ter. We then compare the predictions of the fitting formulae 
with the intrinsic optical depth for SAM model galaxies. We 
find that the proposed fitting formulae are able to provide 



a reasonable fit to Ty in both redshift intervals. The pro- 
posed fitting formulae are analogous to the Guiderdoni & 
Rocca-Volmerange (1987) and De Lucia & Blaizot (2007) 
prescriptions; however, we find different power law depen- 
dencies with respect to those models. 

We compare the GRASIL attenuations with prescriptions 
combining the De Lucia & Blaizot (2007) composite atten- 
uation law with different estimates of Ty. We conclude that 
our new fitting formula is able to reproduce the distribu- 
tion of attenuations and their scaling with galaxy properties, 
such as stellar mass and specific star formation rate, nearly 
as well as the model that uses the actual value of Ty s ex- 
tracted from the RT libraries. The analytic recipes from the 
literature (Guiderdoni & Rocca-Volmerange 1987; De Lucia 
& Blaizot 2007) do reasonably well at low-z but give poor 
results at high-z. 

It is important to keep in mind that the SAM can pro- 
vide many physical properties of galaxies to GRASIL, such 
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Figure 12. Effect of the different prescriptions for dust absorption on the galaxy luminosity function in different bands. Lines as in 
fig. \E\ thin solid line refers to the unextinguished luminosity function. 



as star formation rates, gas density and metallicity, and disc 
and bulge sizes, but it cannot provide many important physi- 
cal quantities of the ISM, like the fraction of gas in molecular 
clouds, the timescale after which massive stars emerge from 
the heavily obscuring molecular clouds or the distribution 
of grain sizes. For the sake of simplicity, we have kept the 
values of these parameters fixed for all model galaxies in our 
libraries. Our choice reflects a parameter combination that 
has been proven adequate for reproducing the properties of 
z < 3 galaxies (Silva et al., 1998; Fontanot et al., 2007). 
Little is known about the redshift evolution of dust grain 
properties: a strong systematic variation of these parame- 
ters may indeed affect the validity of our fitting formulae. 
Observations of dust properties in the host galaxies of high-z 
quasars (see e.g., Maiolino et al., 2004) suggest that the in- 
trinsic properties of dust grains (i.e. composition and dimen- 
sion) are indeed likely to change at z > 4. This is probably 
due to the fact that at such redshifts dust is mainly pro- 



duced by Type II supernovae, whereas at lower redshift the 
contribution from the envelopes of evolved low-mass stars 
becomes dominant. Therefore we conclude that our fitting 
formulae may be limited by this effect to z < 3. Finally, it 
is worth noting that even neglecting any variation of dust 
properties with redshift, we still find a large dispersion in 
the resulting attenuation laws: this result suggests that the 
effect of the geometry and the complexity of star formation 
histories dominate over the properties of the dust grains (see 
e.g. Granato et al. 2000). 

Our results show that the use of a spectro-photometric 
code coupled with a simple attenuation law may introduce 
a systematic shift in the predicted magnitudes. This effect 
has to be taken into account as an uncertainty associated 
with models when comparing their prediction to observa- 
tions. The best way to decrease this systematic error is to 
use a full RT computation. 
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Figure 13. Edge-on relative to face-on attenuation, Ago — Ao, as 
a function of stellar mass. The shaded area shows the GRASIL 
results, the vertical texture the prediction of Prescription TAU- 
DLB 
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